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1. Introduction 

Achieving long time simulations of binary black hole spacetimes remains one of the 
most pressing issues in numerical relativity. Numerical simulations are plagued with 
instabilities, the origin of which has not been isolated among a number of possible 
factors of relevance to stability (see for instance 01213] and references therein). A 
partial list of suspected sources of instabilities includes ill-posed evolution, ill-posed 
constraint propagation and poor choices of boundary conditions, as well as poor choices 
of binary black hole data. The problem of well posedness of the evolution equations 
and constraint propagation in the analytic (as opposed to numerical) sense has been 
widely studied (see the Living Review article by Reula 4 ), with the result that there 
is a great deal of choice of strongly hyperbolic formulations of the Einstein equations 
available for analytic studies (aside from numerical implementation). On the other 
hand, appropriate boundary values for evolution remain to be identified. Since a great 
deal more about the boundary value problem of a strongly hyperbolic set of equations 
is known than for any other kind, currently analytic investigations of appropriate 
boundary conditions tend to use some kind of strongly hyperbolic formulation of the 
Einstein equations. Relevant considerations usually run along the following lines. 

The boundary values of a strongly hyperbolic system of partial differential 
equations split always into two sets those that are determined by the initial 
conditions of the problem and cannot be chosen freely, and those that are entirely 
arbitrary because they are determined by values of the fields outside of the region of 
interest. In practice, the free boundary data must be specified in order for a unique 
solution of the equations to be determined in the region of interest. 

Suppose we have a boundary delimiting a region where we seek a solution to 
the Einstein equations in some strongly hyperbolic formulation. In this case one 
calculates data on the initial slice satisfying the constraints, and then integrates up 
the hyperbolic equations step by step. But, by the previous paragraph, the initial data 
alone are not sufficient to find a unique solution in the future of a bounded region in 
space. Boundary values must be prescribed as well. Suppose that we know exactly 
which variables' boundary values are to be considered independent of the initial data. 
How are we to prescribe them? Since in some sense the solution that we are seeking 
must satisfy all the Einstein equations, and by construction the evolution equations 
are being imposed, then it seems that the appropriate question to ask is: What are 
the boundary data that preserve the constraints? 

This question is addressed by Stewart in j^i within a particular strongly hyperbolic 
formulation of the Einstein equations [71|S|. The constraints themselves, considered 
as real functions of spacetime, propagate according to a strongly hyperbolic system of 
their own. This implies that vanishing values of the constraints at the initial time will 
propagate towards the future along characteristics. Some of the constraints propagate 
towards the boundary and cross out of the region of interest, whereas others propagate 
into the region of interest by crossing the boundary from outside. Clearly the values of 
the incoming constraints at the boundary are arbitrary, and one wants to have them 
vanish. But the vanishing of the constraints cannot be imposed along the boundaries 
in practice. The constraints involve derivatives of the fields across the boundaries, 
not just the values of the fields themselves. Stewart argues that the vanishing of 
the incoming constraints can in fact be used as a boundary condition on the fields 
when the equations are linearized around flat space. If the fields are expressed in 
integral form in terms of Fourier-Laplace transforms, the linearity of the differential 
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equations implies algebraic equations for the transforms of the fields. Additionally, 
the constraints transform into algebraic expressions in terms of the transforms of 
the fields, thus acquiring an algebraic look. Regardless of whether these seemingly 
algebraic conditions are practical boundary conditions for a numerical simulation, the 
argument does not hold up in the non-linear case. 

The idea of imposing the vanishing of the incoming constraints as boundary 
conditions is pursued further in P| , where space derivatives of the fields are eliminated 
in favor of time derivatives in the expression of the incoming constraints in terms of the 
fundamental variables. In this different formulation of the Einstein equations 

is used QHIj restricted to spherical symmetry. 

Clearly not all formulations have the problem of imposing the constraints at 
the boundary. Suppose there is a formulation of the Einstein equations that has 
no incoming constraints, that is: a formulation in which the constraints propagate 
upwards along the boundary. In such a case, the constraints will be satisfied at the 
boundary by virtue of the initial values alone, and there is no need to impose additional 
conditions on the boundary values in order to enforce them. The boundary values of 
the incoming fields in this case must be arbitrary. Such formulations exist at least 
in symmetry-reduced cases, and some boundary problems for those formulations have 
been studied [11]. 

Here we stray away from the general trend of seeking a way to impose the 
constraints along the boundary, in order to propose a seemingly unrelated method to 
write down equations that must hold among the boundary values of many quite generic 
first-order formulations of the Einstein equations. These consist of the vanishing of 
the four components of the projection of the Einstein tensor along the normal to the 
boundary. 

In Section |21 we briefly describe the formulation of the Einstein equations that 
we choose to write down the proposed boundary conditions explicitly in a model 
case. The discussion in Section |21 hinges heavily on the existence of a complete set of 
characteristic fields, which is guaranteed by the strong hyperbolicity of the formulation 
chosen (as opposed to weak hyperbolicity). In Section Owe show that such boundary 
conditions are consistent with constraint propagation, realize the goal of guaranteeing 
the vanishing of the incoming constraints at the boundary and, furthermore, they 
coincide with the boundary conditions proposed in [Hj by the trading of space and 
time derivatives alluded to above. The discussion on Section |31 hinges heavily on 
the existence of a complete set of characteristic constraint fields, namely: on the 
fact that the constraint propagation is strongly hyperbolic (as opposed to weakly 
hyperbolic). In Section 0] we state the generality of the arguments with regards 
to three-dimensional strongly hyperbolic formulations, namely: the aspects of the 
argument that apply to any strongly hyperbolic formulation of the Einstein equations 
irrespective of particularities such as the specific characteristic fields and speeds. 
Formulation-dependent features such as any explicit form of the boundary conditions 
arising from the projection of the Einstein equations perpendicularly to a boundary are 
intentionally excluded from Section^jwith the purpose of exposing the true reach and 
relevance of the argument across the lines delimiting alternative strongly hyperbolic 
formulations. Concluding remarks appear in Sectional 
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2. Boundary conditions for the Einstein-ChristofFel formulation with 
spherical symmetry 



As usual, in a foliation of spacetinie by level surfaces of a time function t, and using t 
as a coordinate, the spacetime metric can be viewed as a time dependent Ricmannian 
metric evolving in time according to a freely specifiable rate and labeling of the 3- 
space. From now on, we use the term metric with no qualifier to refer to the evolving 
Riemannian metric of the slices. In this section we use the Einstein-Christoffel (EC) 
formulation jlU|. restricted to spherical symmetry. The full 3-D equations and their 
structure of characteristics can be found in ^^li as well as the restriction to spherical 
symmetry. We borrow the notation directly from TT . With the restriction to spherical 
symmetry, the spacetime metric has the form 

ds^ = ^N^dt^ + g„{dr + dtf + grr^idO^ + sin^ Bdcf)^) (1) 

where grr,9T^ and N are functions only of r and t. The shift /?'" and the densitized 
lapse 

" (2) 



9Ty/grr 

are considered arbitrarily given. The EC formulation of the Einstein equa- 
tions restricted to spherical symmetry consists of six fundamental variables 
(grr, gr, Krr, Kt, frrr, frr) which cvolvc according to six evolution equations and 
whose initial data satisfies four constraints I12j . The evolution equations are 



dtgrr — drgrr — 

dtgr - P^'drgr = 
N 

dtKrr — P^drKrr H drfrrr = N 

grr 



2NK^ 



2grrdrP'', 

13'' 



2NKt + 2—gT, 
r 



(3 a) 
(36) 



2/ rr [ f rr ~\~ 



4.f, 



rT 



9T 



Kr 



Kt 
I 

' gr 



- ^ — 6 ( ) — dt In a 



6 



frrr P fr 



OrJrT 

grr 



H OrKr 

grr 



r 



+ ATTN{Tg.r 
N {ktK - 
2/3'- 



dr In a 
2/, 



{drXnaf 



1 

l2 



rT 



frT rs , ~ 

Or m a 



grrgr gr 



-Kj 



(3c) 



(3d) 



Kt fj, 



gT 



Agrr— 3^ - r 



gT 



gT 



(?r In a 

r 



9^ In a 



-t- ifrrrdrP'' + ffrr^,^/?'' + IGnNJrC, 



N 

dtfrT - P^'drfrT H OrKr 



= NKt 



r.frT _ rr 
^ J r 

gT 



9r In ( 



drP"- + ^ ) frT. 



(3e) 



5 



The constraints are 



9rr9T '^r^gt grrQT \r 2gT J gr \ 'igrj 

C .Ml + ^ + /:lf^; + ^U4... = 0, (46) 



5t rgr gr \ gr J 

r = rl J- ^9rrfrT 9 f r _ n (A^ 

CrT = drgr + ^ - 2frT = 0. (4d) 

r 

Here C and Cr are the scalar and vector constraints, respectively, whereas the vanishing 
of Crrr and CrT defines the variables frrr and /,t- We are here interested only in the 
vacuum equations, so we set all the source terms to zero (T — Srr ~ Jr = P = 0), and 
we do not need to consider matter equations. 

In essence, the evolution system is strongly hyperbolic because one can find a 
complete set of linearly independent combinations of the fundamental variables such 
that the principal part of the equations decouples. Such combinations are referred 
to as the characteristic variables. In a linear homogeneous system, the characteristic 
variables propagate their initial values exactly along the characteristic lines, and may 
be interpreted as waves propagating with the characteristic speed. In our case |12| . 
the characteristic fields and their characteristic speeds are 

= Srr (fc = PI (5 a) 

= ST («c = n (56) 
C/± = Krr ± 4^ {v, = (3''T agr) (5c) 



U^ = Kt± ^ {vc = T agr) (5d) 

Y grr 

Here we have changed the signs of the characteristic speeds relative to JT? in order 
to agree with the sign conventions in 0. The characteristic fields give insights into 
how much data one can or must provide in addition to the initial data, depending 
on the region of spacetime where the solution is sought. If the region of interest has 
a boundary, in general one may need to prescribe additional data on this boundary 
in order to obtain a unique solution within that region. The additional data are the 
values of the characteristic fields that enter the region from outside. If such values 
are not prescribed, then there are many solutions for the same initial data in the 
region of interest. For each boundary value prescribed arbitrarily - as long as it is 
consistent with the initial data at the intersection of the boundary with the initial 
slice -, there is a unique solution. No boundary value may be freely prescribed for 
the characteristic fields that cross the boundary from the inside, or that run along the 
boundary, since such values are determined by the initial data already, in principle. 
In practice, since the equations are not linear nor homogeneous, the actual values of 
the outgoing characteristic fields at the boundary are not known in advance. Most 
of the time the values of some or even all the fields at the boundary, regardless of 
formulation type, are prescribed by means of some type of wave condition, as in 
and, more elaborately, in |14j . 

In our problem, we restrict attention now to a region inside a fixed value of r. 
There is a set of incoming characteristic fields, which are the ones that have positive 
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characteristic speeds in our convention, which we wiU be able to single out once we 
choose the values of the densitized lapse and shift. Suppose initial data are prescribed 
with vanishing values of the constraints. One can now prescribe the values of the 
incoming fields independently of the initial values. The question is: are such values 
free? It is worth investigating the possibility that the Einstein equations themselves 
may restrict the boundary values. 

Consider the following argument. The boundary is a surface of constant radius 
in spacetime. The unit normal vector to the boundary is well defined. We denote it 
by e" = (e*, e'', e^, e"*^) and it can easily be computed from the gradient vector r,a = 
(0, 1,0,0) by e° ='*g°V,b /\/V%^ =^ 9"'' / VY^ , where Y'' is the contravariant 
spacetime metric. In our case, we have 



or 

e* = , ^ , =, (6a) 



agT^/grr 

e» = = 0. (6c) 

The projector pat on the boundary surface can be found from the spacetime metric 
and the unit normal to the surface via 

Pab = %ab - eaCb (7) 

with latin indices a,b,c . . . raised and lowered with the spacetime metric. In particular, 
we are interested in the projector with mixed indices 

p\ = S\ - e^et (8) 

We saturate the indices of the Einstein equations Gat = in two different ways. 
One equation is obtained by contracting twice with the normal to the boundary 
surface, e°. We have 

Gabe^'e'' = Gttie'f + 2Gt,e'^e* + Grrie^^ = (9) 

Another set of equations is obtained by contracting one idex with 6°" and the 
other index with p^c, namely: Gabs'^P^c = 0. But because of the spherical symmetry, 
this contraction is identically zero for c — 0,(j). The remaining equations are 

GabeVt = Gtte' + Gt.e'' = 0, (10) 

Gabe'^pK = {Gtte* + Gtre^pK + (Gt,e* + GrrC^pW = 0. (11) 

Clearly, if H10|l and JHl are satisfied, then Ijllll is an identity. Thus there are 
only two independent equations, and we pick (|l()|l and Q. These equations are 
extremely valuable for the boundary value problem because they contain no second- 
order derivatives with respect to r irrespectively of the choice of or a. In the case 
of a first-order formulation, like ours, such equations are said to be interior to the 
boundary surface - not to confuse with the interior of our region! In the following, 
we restrict ourselves to /?'" = and 5=1 just for the sake of argument. Non trivial 
choices of densitized lapse and shift only make the argument more complicated without 
introducing any real obstacle. 

Explicitly, for the case of = and 5 = 1, we have 

Gabe'^Cb = = 5 \ 5- (2r^grr^^^kT - Sr^grrfrT^ " "^rgrgrrfrT 

r^gr grr V 
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+ 2r'^gTfrTfrrr - 9rr^ gx + r^grrKr^j , (12a) 

b 2 / • 

Gabe^P t = = :rpr rg„frT + rgT\/g^frTKrr - ^/g^Kr {rgrrfrT 

rgxgrr'^ ^ 

+2gTgrr - rgrfrrr) ) • (12&) 
We can put these equations entirely in terms of the characteristic fields, using 

At =1{U} + Ut), (13a) 
frT =^{U+- U^) , (13&) 



Krr = l{U+ + U-), (13c) 
frrr = ^ {U+ - U') . (13d) 

and the fact that the metric components grr and gx are characteristic fields themselves. 
For the derivatives Kt and frT we take derivatives of p3a|l and (|136|) and use H3a|l 
- written in terms of the characteristic fields- to substitute grr which will appear in 
frT- Explicitly, for = and 5 = 1 we have 

Kt = 1 {u} + Ut) (14a) 

9t 



f,r ^^[U+- U^) ^ (C/+ + U-) (C/+ - U^) (Ub) 

On account of p4Qll4b|l . one can see by inspection that when both equations (|12at 
I126|) are written out entirely in terms of characteristic fields, they involve the time 
derivatives only of U^, but not of U^. Furthermore, appear linearly in both 
equations. We can algebraically solve the two equations for both tl^ as independent 
variables in terms of the rest. We readily find: 

1 / 2 

U^^ r, 2 /= [9rrgT + r'^grr (Uj^) - 2r'^ grrUj^ U:^ 

- 'ir^gTU^+r^gTUTU+- r^gTUj^Ur) (15) 



and 



1 / 2 

U^= „ n \9rrgT + r'^grr (C^T ) — 2r'^grrU^Uj 



'^-2r^^rr 

+ ^r^gTU+- r^gTU+U++r^gTU+Ur^ (16) 

Both equations H15|) and (|16() are necessary because they are linearly independent. The 
question is: how are we to interpret them within the framework of the initial-boundary 
value problem? 

To start with, the fields grr and may be considered as known sources, because 
they travel upwards along the boundary. If the problem was linear and homogeneous, 
then the values of grr and gT at any point on the boundary would be exactly their 
initial values. Because the problem is non-linear, the values at any point on the 
boundary must be integrated, but can still be regarded as sources. A similar thinking 
may be used for the outgoing characteristic field ?7+, with the complication that C/+ 
travels towards the boundary but not directly upwards, and thus the computation of 
its values at the boundary will be more involved, but conceptually no different. 
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Additionally, we may take the point of view that the incoming characteristic field 
U~ is arbitrary, since there are only two equations to satisfy at the boundary and 
U~ does not appear in any of them. From this perspective, the value of this field at 
the boundary is truly free, and must be considered as one true degree of freedom of 
the boundary value problem for the Einstein equations, in addition to the degrees of 
freedom contained in the initial data. 

Continuing within such an interpretation, Eq. Ijl5() provides exact boundary values 
for the incoming characteristic field in terms of the free incoming field U~ and 
the fields that either run along the boundary or cross from the inside, all of which are 
determined by the initial data in principle, as argued before. 

There remains Eq. (|16|l . From our point of view, this equation predicts boundary 
values for . There is no doubt that it must hold, in order for the solution that 
we are seeking to satisfy all the Einstein equations at all points on the boundary 
surface. But is, in principle, determined by the initial data by propagation along 
characteristics, as explained before. There seem to be two possibilities here: either 
the outgoing field propagated from the initial data satisfies Eq. ((TB|l . or it does 
not. If the boundary values of do not satisfy Eq. then the initial-boundary 
value problem is inconsistent and can not yield a solution to the Einstein equations 
in the region of interest. 

Alternatively, if the boundary values of propagated along characteristics do 
satisfy Eq. (|16|l . then the initial-boundary value problem is consistent. In this case, 
Eq. (|16|l can be used to prescribe the values of with more accuracy and reliability 
than the propagation along characteristics discussed above because it is an ordinary 
differential equation, and should be preferred. 

In the following Section we show that the initial-boundary value problem is 
consistent and we should be using Eq. (|16|l to prescribe the boundary values of 

instead of propagating them by characteristics from the initial data, because the 
propagation by characteristics is inaccurate in itself -even more so in the case that 
the characteristic speeds depend on the fields themselves and must be calculated at 
every step. One can anticipate that this is the case based on some intuition. Why 
should the boundary values of the outgoing fields be "constrained" by an equation such 
as Eq. ((TB|l ? Because, in principle, the boundary values of the outgoing fields must 
reflect in some way the constraints on the intial data that determine them. The initial 
data are related; their inter-relationships must be propagated along characteristics. 
Eq. (|16(l may well be a shortcut out of constraint propagation. 



3. Consistency with constraint propagation 

If thought of as functions which could take any real values, the initial constraints 
C,Cr,Crrr,CrT cvolve in time according to another strongly hyperbolic evolution 
system, with the following characteristic fields and speeds 

- agr) (17a) 
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In the case of /J'' = and a = 1, C2 is the only incoming characteristic field. If 
one wishes for all four functions to take the value zero in the region interior to some 
fixed radius, then three of them will be zero by choice of initial data, but C2 must 
be set to zero at the boundary in order to propagate inwardly from there. But one 
simply may set the value of C2 to zero because one does not evolve the constraint 
propagation system. It is one thing to impose the value on the function C2, but 
a different thing to expect C2 to be vanishing at the boundary as a function of the 
fundamental variables of evolution. In particular, by inspection one can readily see 
that C2 contains the derivative of with respect to r. Therefore, although it is 
natural and desirable to have C2 = on the boundary, it is an impossible task to have 
C2 = as a boundary condition for the fundamental fields of evolution. C2 will have 
to vanish for the evolution to produce a solution of the Einstein equations, and thus 
the vanishing of C2 and all the other constraints must be checked for consistency at the 
boundary, after the solution has been found and the r— derivatives can be evaluated; 
but it is not practical as a boundary condition. 

A strategy that has been proposed [5] to circumvent this obstacle is to use the 
evolution equations to turn the r— derivatives of the characteristic fields that appear in 
the expression of C2 into derivatives. The procedure must work in this symmetry- 
reduced case because, by construction, all characteristic fields have time derivatives 
proportional to their r— derivatives, up to terms of zeroth order. Therefore, C2 = can 
be turned into an evolution equation for that will be restricted to the boundary. 
Following we do that next, maintaining the restrictions (i^ — Q and 5 = 1. We 
have explicitly 

'^2 = " \ ^ [2r^gTgrrdrKT " grOrr^^^ Or frT + grgrr^'"^ 

gT grr \ 

- Argrgrr^^^ frT - 7r'^grr^^^.frT + 2r'^ gr y/g^ frT frrr + 2r^ gr grr^^^ Kt Krr 

+ r2g,//2ifj,2 ^ ^rgtgl^Kr - 2r^gTgrrfrTKrr - 2r^glfrTKT^ (18) 

Solving for drfrT from the evolution equation (|3rf|l . and for drKx from the evolution 
equation ( |3^ , substituting into IjlHI) . and replacing all appearances of frT, frrr, Kt 
and Krr in terms of the characteristic fields using H13all3cj|) , Eq. (|18|l turns into a new 
expression, which is not an initial constraint anymore, so we use quotes to make this 
point explicit: 

"^2" = - ^ 2 2 5/2 (^r^9rr^U^ - r'^grr'^'HU^f ^ grgrr^'^^U^ U+ 
2r'^gT grr ' ^ 

+r^gTgrr'"^U^U- - grgrr^/^ + 2r^grr^/^U+U^ + ^rgrgrr'^U^) (19) 

Setting "C2" = yields a boundary condition for . The point is that C2 = 
on the boundary is not the same as "C2" = 0, but they are the same along the 
boundary if evaluated on fundamental fields that satisfy the evolution equations. In 
other words: they diff'er by a linear combination of the evolution equations. But more 
interestingly, "C2" = is exactly the same equation as (|15f) . as can be verified by 
inspection. Therefore the method of "trading" space derivatives for time derivatives 
advocated in W is equivalent to solving one of the Einstein equations that are interior 
to the boundary, there being a two-dimensional set of such equations. This could 
be expected because the projections of the Einstein equations along the direction 
normal to a surface of fixed radius are the linear combinations that have no second- 
order r— derivatives, and clearly the "trading" is equivalent to linearly combining the 
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Einstein equations. Now, because Eq. (|16|l is equivalent to C2 = on the boundary, 
it guarantees that C2 wiU vanish in the interior by propagation along characteristics. 
So Eq. (|15|) is not only consistent with the boundary value problem, but is necessary 
as well. 

We can also use the "trading" of space derivatives for time derivatives to write 
an equation "Ci" = from Ci = along the boundary. This is exactly Eq. 116|l . the 
meaning of which is explained in the previous section, and we conclude that Eq. (|16|l 
and Ci are equivalent along the boundary if the evolution equations are satisfied. But 
how much of this argument proves that the constraint has been "propagated" and 
where does the propagation start? Clearly Ci propagates out from the initial data 
towards the boundary, therefore Ci will vanish along the boundary as a consequence 
of its vanishing initially. Since the evolution will be satisfied by construction, then 
"Ci" will vanish as a consequence of Ci vanishing initially as well. Seemingly, thus, 
the boundary values of the fields satisfying Eq. (|15|) must be consistent with the 
initial values satisfying Ci = 0. For this reason, and because Eq. (|15|l is an ordinary 
differential equation, it is not only consistent but advisable to use Eq. (|16|l to obtain the 
values of instead of propagating them from the initial values along characteristics. 
At any rate, Eq. Ijlfifl should be checked for consistency regardless of method. 

Two objections may come to mind on a surface glance to Eq. (|16|l . which we 
want to anticipate. It might conceivably be argued that if the values of C/^ at the 
boundary are obtained from the initial values by propagation along characteristics, 
then Eq. (fT?)|l could also be thought of as an equation for the incoming field U~ 
with given, instead of prescribing with a free U~ . But this is a flawed line of 
thought, for the boundary values of propagated along characteristics are consistent 
with Eq. H16|l . and substituting into Eq. Hlt)|) must yield an identity, leaving no 
equation to solve for U~ , which would thus remain arbitrary. 

Secondly, how is it that Eq. H16(l involving a free boundary field U~ can be 
consistent with propagation of along characteristics, given that C/^ is a field which 
the initial data know nothing about? The initial data do not know about U~ , but 
the characteristics do, because the problem is non linear. Thus the propagation of 
along characteristics will necessarily yield boundary values of Uj^ that know about 
earlier boundary values of incoming fields. It is fine for Eq. H16|) to involve U~ , as 
long as it does not involve il~ . 

For additional support of our scheme to use Eq. l(T5|l and Eq. fTH|l to prescribe U,^ 
and , respectively, assuming that U~ is given freely, we can linearize the equations 
around flat space and verify that the scheme is consistent, since the boundary value 
problem of linear hyperbolic equations is clear and enlightening. We have 



Qrr 


^ 1+grr 


(20a) 


9T 


= 1 + ?T 


(20&) 


Krr 


— 


(20c) 


Kt 


= Kt 


(20rf) 


I rrr 


4 

— frrr 
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(20e) 


frT 


= -+frT 

r 


(20/) 



where all hatted quantities are small. We will keep only linear terms in such quantities. 
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The evolution equations become 

grr = - ^Krr (21a) 

= -2Kt (21 &) 

c, 7 42_ 48_ lOj- 44- 

Krr + Orfrrr — n9rr H TffT H frrr JrT, (21c) 

i^T + Sr/rT = + - -frT (21(i) 

+ a^i?^^ = - yKrr + ^Kt (21e) 

tr + drKr = - -Kt {2lf) 
r 

Thus the characteristic fields are , Qt and 

U^ = Kt± frT, (22a) 

= Krr ± frrr- {22b) 

Our boundary conditions, Eqs. (|15|l and l|16l) . Unearize to 

= ^9rr - ^gr + + ic/+ , (23) 

and 

9 ^ 3 ^ 1 

respectively. Thus (|23|l can be used to prescribe values for Uj^ if U^T is given arbitrarily, 
whereas can be used to calculate values for Ui^ irrespective of U~ , as anticipated. 



Ut = - T^gr - -U+ ~ -U^ , (24) 



4. Boundary conditions for generic three-dimensional strongly hyperbolic 
formulations 

How much of the argument in spherical symmetry actually depends on the symmetry 
assumption? Not a great deal. Suppose we have a strongly hyperbolic formulation 
of the Einstein equations in terms of 6+6+18 variables representing the three- metric 
and all its first derivatives. The evolution equations then look like 

ii = A'diu + b (25) 

where u is the 30-dimensional vector of all the fundamental variables, and there are 
4+18 constraints on the initial data: 

C =0 (26a) 

C, ^ (266) 

C^Jk = (26c) 

where C and Ci are the scalar and vector constraint, respectively, and Cijk are the 
constraints necessary to reduce the equations from second to first order in space (they 
define the 18 first order variables). Wherever the set of evolution equations and the 
constraints are satisfied, the ten Einstein equations Gab = for the ten components 
of the spacetime metric gab are satisfied, equivalently. 

Suppose the evolution equations are strongly hyperbolic, and that they imply, 
in the usual manner jl5| . a second (or subsidiary) system of equations for the constraint 
functions C,Ci and Cijk- Suppose that this second system of equations, for the 
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constraints, is also strongly hyperbolic. Assume we have identified the characteristic 
variables of both strongly hyperbolic systems. 

Suppose now that we seek a solution in the region interior to some fixed value 
of the coordinate (any spacelike coordinate). The unit normal to the boundary is 
then 

= ^ = ^ (27) 



and the projector on the boundary surface is = ^1 + e^e^ with Cf, — 51/ \/ g^^- We 
can write down Gabe°'e^ — and Gabe'^p^c — 0. In general these will be four equations 
with no second derivatives with respect to x^, as can be proven by direct calculation. 
Contracting the Einstein tensor with e" yields a vector equation with two pieces: 

Gate'' = Rate'' - -^R- (28) 



We next show that the components 6 7^ 1 of Rate'^ have no second derivatives of 
the metric with respect to a;^, and that the Ricci scalar term cancels out the second 
a:^— derivatives that appear in the first term for b = 1. We use the following well- 
known expression of the Ricci tensor |16 | in which the second-order derivatives of the 
metric are explicit: 

Rab — 2-9'^'' {9cb,ad + gad,cb — Qcd^ab — gab,cd) + g'^'^ (Xld^ edb ^ ^ab^ecd) ■ (29) 

Contracting with yields 

D „a / la Ic 11 cd \ , io(\\ 

^'^^e = —7= [g g gacib - g g gcd.ib) + ■■■ (30) 

2V5 

where . . . represents terms that have no second derivatives with respect to x^ . Thus 
Rabe'^ has no second a: ^—derivatives except for 6=1. The second part of the Einstein 
tensor is nonvanishing only for 6=1. Therefore the only component of Gabe" that 
might contain a second x^— derivative is 6 = 1. Now R is exactly 

R = 5'"5''5ac,ii - g''g''gcd,ii + ■■■ (3i) 

Therefore, the terms with second a; ^—derivatives of the Ricci scalar part of the Einstein 
tensor cancel exactly with those in the Ricci tensor part, with the consequence that 
Gabc"" has no second a; ^—derivatives of the metric for any value of 6. Contracting 
with e'' or with p^c will not add new second x^— derivatives, but provides perhaps a 
convenient splitting of the four equations. The four equations 

= Gabe'' = (32) 

or equivalently 

G„6e'^e'' = and GaheV ^ = (33) 

must be satisfied by any solution to the Einstein equations and are interior to the 
boundary surface, which makes them ideal boundary conditions. There are four 
equations and thirty variables, but of the thirty variables many are outgoing or run 
along the boundary. In fact, one can always write the evolution system so that the 
metric components propagate with vanishing characteristic speeds, regardless of the 
choice of lapse and shift. So we have at most (30— 6) /2 relevant incoming characteristic 
fields. Clearly, some of the incoming fields will be left arbitrary at the boundary - 
perhaps more than 8, since in our example in spherical symmetry we find that one of 
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the boundary conditions applies to an outgoing field instead of an incoming one. The 
particulars of which characteristic fields will be affected by the boundary equations 
will depend entirely on the details of the first-order strongly hyperbolic formulation 
of choice. As an illustration, the direct generalization of our argument to the case 
of the Einstein-Christoffel formulation in three-dimensions without the restriction of 
spherical symmetry will be reported elsewhere. 

Conceptually, the scheme runs as follows. Choose initial data that satisfy l@J. 
Choose boundary data that satisfy 1)32(1 . The vanishing of the initial constraints 
guarantees that the outgoing and static constraints will remain vanishing in the 
region of interest. The boundary equations (|32|l are equivalent to the vanishing of the 
incoming constraints at the boundary on account of their being related to the incoming 
constraints by linear combinations with the Einstein equations that are satisfied at 
the boundary - the evolution is satisfied by construction and the outgoing constraints 
are satisfied by propagation along characteristics. Therefore the boundary equations 
H32|l guarantee that the incoming constraints are vanishing at the boundary, which in 
turn guarantees that they will remain vanishing in the region of interest. 

Although the preceding calculation assumes that the normal to the surface of 
fixed value of is spacelike, clearly Eqs. ((32(1 will be interior to the surface of fixed 
value of even if the normal is timelike, in which case the proof that the equations 
contain no second x^— derivatives runs the same if we substitute with y^jg^. 

This remark applies to interior boundaries that lie within the event horizon of a black 
hole spacetime, and may be relevant to numerical simulations of binary black holes. 

On the other hand, in this Section we have restricted the discussion to three- 
dimensional boundaries taken one at a time. As is the case with any three-dimensional 
boundary- value problem, the intersection of two boundaries requires special attention. 
In this respect, consistency issues between the boundary equations that may arise at 
the intersection of two boundary surfaces remain to be studied, and their resolution 
may depend on the particulars of the formulation of the initial value problem at hand. 
The interested reader is referred to [17j where a particular representation of the 
projection of the Einstein equations on the boundary is implemented and corners and 
edges are treated, for the case of an ADM-like formulation of the initial value problem. 

5. Concluding remarks and outlook 

Even though we have, for the most part, developed our arguments explicitly in the 
case of spherical symmetry for the EC formulation, the fundamental relevance of 
the argument to the initial-boundary value problem of the Einstein equations does 
not depend on the symmetry restrictions nor on the particulars of the hyperbolic 
formulation, as we show in Section 0] 

We argue that given any boundary at a fixed value of a coordinate of an initial- 
boundary value problem for the Einstein equations, the vanishing of the components 
of the projection of the Einstein tensor Gab along the normal to the boundary, 
namely Gate^ = = 0, constitute necessary and consistent boundary conditions, for 
essential reasons, as follows. 

First, the components of the projection of the Einstein tensor along the normal 
to such a boundary, Ba, contain no second-order derivatives sticking out of the 
boundary. If the initial value problem is stated in first order form, such as any of 
the hyperbolic formulations available, their vanishing becomes differential equations 
for the boundary values of the fundamental fields excluding their derivatives across 
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the boundary. Moreover, if the initial value problem is stated in second-order form 
in the space coordinates, such as ADM form ^l^j or conformal form ^11201, they 
constitute a type of mixed Neumann-Dirichlet conditions, the appropriateness of which 
is worth investigating in full depth. In this regard, some representations of Sa = 
as boundary conditions have been used along with an ADM type formulation of the 
linearized Einstein equations in \17\ to investigate the consistency of pure Neumann 
or Dirichlet conditions with the evolution. 

Second, in the case of strongly hyperbolic formulations that propagate the 
constraints in a strongly hyperbolic fashion, the vanishing of all four boundary 
equations Ba guarantees the vanishing of the incoming constraints in the region of 
interest and is consistent with the propagation of the outgoing constraints along 
characteristics. In practice, this means that one can use ;Ba = to prescribe 
as many boundary values as possible, without worrying about inconsistencies with 
the initial values. This necessary analytic consistency should not be mistaken for 
numerical consistency, however. In fact, it is not known at this stage whether this 
analytic consistency is stable under small perturbations of either the initial data or the 
boundary data, a point that is critical to the numerical implementation. Additionally, 
how to proceed in order to investigate the stability of an initial-boundary-value 
problem with constraints remains, to our knowledge, an open problem not addressed 
in the standard reference literature of strongly hyperbolic systems such as 

In the case of formulations that do not guarantee stable constraint propagation, 
= on the boundary still appeals to us as the next best thing. In fact, in such cases 
there are hardly any handles on boundary conditions. Notice that imposing Ba — on 
the boundary does not involve any manipulation of the choice of evolution equations; it 
does not affect the evolution equations themselves nor the system of propagation of the 
constraints. The converse is not true: any manipulation of the evolution equations (by 
adding linear combinations of the constraints in the usual manner) necessarily affects 
both the propagation of the constraints and the role of i3a = on the boundary. 
The interested reader is referred to for a series of numerical studies of the effect 
of the manipulation of the evolution equations on constraint propagation exclusively, 
without imposing Ba — 0. 

Several issues remain open at this time. First and foremost, there remains the 
issue of whether the boundary equations Ba ~ are consistent with a well-posed 
initial-boundary-value problem in the analytic sense, that is: in the sense that the 
solutions at a later time are continuous functions of the initial data and the boundary 
free data. The answer to this question will vary among the different strongly hyperbolic 
formulations available. The answer may be relevant to numerical relativity because, 
even though well-posedness of the initial-boundary value problem does not guarantee a 
stable numerical implementation, it is considered as a necessary or at least desirable 
feature of the continuous equations being implemented. Second, any issues that may 
arise exclusively in connection with the numerical implementation of Ba = -but 
that are otherwise irrelevant to the analytic initial-boundary value problem of the 
Einstein equations- are also open at this time and are worth pursuing, their details 
being strongly dependent on the particulars of the formulation at hand. From the 
numerical point of view, for instance, the presence of undifferentiated terms, -which 
is not relevant to our current argument and does not affect the well-posedness of 
first-order problems- is critical to numerical stability [S]. In this respect, the effect 
of the undifferentiated terms of the boundary equations ,Ba = on any issues of 
stability remains to be determined, being, almost certainly, critically dependent on 
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the particulars of the formulation at hand. Additional open issues at this time include 
checking for the consistency of = with other boundary conditions already in use 
in numerical simulations, and, eventually, how the use of Ba ~ on the boundary may 
affect the run time of numerical simulations. 
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